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Abstract 

We present evidence that the concurrent existence of two populations of particles with different 
effective diameters is not a prerequisite for the occurrence of anomalous phase behaviors in systems 
of particles interacting through spherically-symmetric unbounded potentials. Our results show 
that an extremely weak softening of the interparticle repulsion, yielding a single nearest-neighbor 
separation, is able to originate a wide spectrum of unconventional features including reentrant 
melting, solid polymorphism, as well as thermodynamic, dynamic, and structural anomalies. These 
findings extend the possibility of anomalous phase behavior to a class of systems much broader 
than currently assumed. 
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I. INTRODUCTION 



The phase behavior of one- component substances is termed "anomalous" whenever it dif- 
fers from that of prototypical (i.e., argon-like) liquids. Anomalous phase behavior includes 
polymorphism in the liquid and solid phases, reentrant melting (i.e., melting by compres- 
sion at constant temperature), and a host of other thermodynamic, dynamic, and structural 
anomalies. An important class of systems displaying such features is that of network-forming 
fluids, i.e., fluids that form orientation-specific, intermolecular bonds which^are strong rela- 
tive to London forces. The most well-known of these substances is water 



cases are silicon, phosphorous, and silica 



% |5[ , while other 
Anomalous behaviors have been observed 



also in systems where effective interatomic forces can reasonably be assumed non-directional. 
For example, alkali and alkali-earth metals at high pressures exhibit one or more regions of 
reentrant melting Jq]. Liquid or amorphous polymorphism have been found in other systems 



10| , triphenyl phosphite , 



with non-directional interactions, such as molten AI2O3-Y2O3 
and some metallic glasses [12 ]. 

The quest for understanding the basic mechanisms of anomalous phase behaviors has 
given impulse to the study of simple isotropic potentials that are able to display the same 
behaviors. These systems can be appropriate models for the generic thermodynamic behavior 
of pure metals, metallic mixtures, electrolytes, and colloids. It has been found that unusual 
behaviors may arise in systems of spherical particles (simple fluids) where the unbounded 
repulsive core is "softened" through the addition of a finite repulsion at intermediate dis- 
tances, so as to generate two distinct length scales in the system: a "hard" one, related to 
the inner core, and a "soft" one, associated with the soft, penetrable, component of the re- 
pulsion 13M26I] . Due to this feature, such core-softened (CS) fluids are characterized by two 



competing, expanded and compact, local arrangements of particles. This property mimics 
the behavior of the more complex network-forming fluids, where loose and compact local 
structures arise from the continuous formation and disruption of the dynamic network orig- 
inated by orientational bonds. So far, the existence of two competing local structures has 
been regarded as essential for the occurrence of anomalous phase behaviors in simple fluids. 
In this work we challenge this belief and provide evidence that these behaviors may also oc- 
cur in systems with a single local particle arrangement, i.e., under much weaker conditions 
than admitted so far. 
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II. THEORY AND SIMULATION 



For a spherically-symmetric unbounded potential u(r) the analytic condition for core 



softening 



27] is that, in a range of interparticle distances r within the core, the product rf(i 



with f(r) = — dw/dr, decreases as r gets smaller. Systems satisfying this core-softening 
condition are characterized by two different repulsive length scales and, for pressures where 
the two scales are both effective, loose and compact local arrangements compete with each 
other and the fluid behaves like a mixture of two particle species with different radii ( "two- 



state" fluid 



28l-l30l|). A typical property of such fluids is the existence of a reentrant melting 



line, a feature that is strictly related to other anomalous behaviors 29|. We note that 

tened two-body 



the core-softening condition is not applicable to lattice models with a so: 
repulsion, which can nonetheless show a reentrant-melting behavior 3ll. |32|. 

In order to assess whether or not the existence of two distinct local structures is a neces- 
sary requirement for reentrant melting, we ask the following question: how much soft should 
repulsion be in order that reentrant melting may occur in simple fluids? To attempt an 
answer, we have developed a convenient tool for the analysis of trends in the melting-line 
topology, based on linear-elasticity theory [331 ] . We focus on the face-centered cubic (fee) and 
the body-centered cubic (bec) phases, i.e., those relevant for purely- repulsive potentials at 
low to moderate pressure. According to elasticity theory, the enthalpy of a cubic-symmetry 
solid under homogeneous strain is 

1 
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Vn 
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where Hq is the enthalpy of the undeformed crystal with volume Vq, e a p are the strain-tensor 
components, and the A's are Lame coefficients. At zero temperature, these coefficients are 
given in terms of the interparticle potential |34j, being ultimately a function of pressure. 
Then, we employ Eq. ([I]) in a field theory where the field variables are the components of 
the displacement vector d(R) at each lattice site. In practice, any given realization of the 
displacement field is assigned the weight exp{— H/iksT)}, with the A's fixed at their T = 
values. In Fourier space, this eventually gives the mean square displacement as a Gaussian 
integral. We find 

(2) 



d2 >=Wos(^ Tr|M " (q)l 
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where the sum in (jSJ) is over the Debye sphere (DS) and M x (q) is the inverse of 
M Q/3 (q) = [A 3 g 2 + (Ai - A 2 - 2\ 3 )q 2 a }8 a p + (A 2 + X 3 )q a q^ . 



(3) 



Derivation of Eq.@ assumes that the quadratic form in ([I]) is positive definite, meaning that 
the stability conditions Ai + 2A 2 > 0, Ai — A 2 > 0, and A3 > must hold. At those pressures 
where any of these inequalities fails to be satisfied, the solid is mechanically unstable. 

We obtain a closed- form expression for the melting temperature T m (P) by using Lin- 
demann's criterion, i.e., by equating the r.h.s. of Eq.(E) for T = T m to (^l^nn) 2 , where 
5l ~ 0.15 for a fee solid and 5l ~ 0.18 for a bec solid [351 ]. and tnn is the nearest-neighbor 
(NN) distance. We have checked in known cases (Lennard- Jones, exp-6, and Gaussian po- 
tentials) that the present criterion correctly predicts the overall trend of the fee and bec 
melting lines, although T m is systematically overestimated, which is not surprising in view 
of the assumption of temperature- independent A's. 

We now apply the above outlined melting criterion to two different types of repulsive 
interactions with tunable softness. We first consider the modified inverse-power (MIP) 
potential WMip( r ) = w (r)[l + Aexp{ — 10(1 — r/cr) 2 }} with u (r) = e(r/o")~ 6 , e and a being 
energy and length units. As \A\ increases, core repulsion becomes softer and the core- 
softening condition is satisfied for A < —0.60 and A > 28. Our criterion predicts reentrant 
melting for A < —0.35 (Fig. 1) and A > 1 (not shown). We perform similar calculations for 
a different interaction law, UYx( r ) = eexp ja(l — r/cr) — 6(1 — r/a) 2 ln(r/eQ j, introduced 
long ago by Yoshida and Kamakura (YK) 36] . This repulsion becomes softer with decreasing 
a and satisfies the core-softening condition for a ^ 2.3. The theoretical melting line shows 
a reentrant portion for a 5.5 (see Fig. 2). 

The above results intriguingly suggest that reentrant melting (and presumably also other 
anomalous behaviors) can occur for continuous potentials that do not satisfy the core- 
softening condition. To assess this possibility, we investigate through numerical simulation 
the equilibrium properties of a system of particles interacting through wykM for a = 3.3. 
This potential does not satisfy the core-softening condition and monotonously increases, for 
decreasing r, together with its first and second derivatives (see Fig. 3 inset). To carry on the 
simulation study, we assume that the same crystals that are stable at T = also give the 
underlying lattice structure for the stable solid phases at T > 0. Through a total-energy 



calculation similar to that in 



37|, we find that at T = the sequence of stable crystals for 
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FIG. 1: (Color online). Theoretical melting lines (see text) for UMip(r). Pressure P and temper- 
ature T are in reduced, e/<r 3 and e/ks, units. Melting lines are obtained by taking, for each P, 
the higher T m value between fee (solid line) and bec (dashed line). Four cases are shown: A = 
(blue), —0.2 (cyan), —0.35 (magenta), —0.55 (red). Repulsion softness grows in the arrow direction. 
Cusps at estimated fluid-fcc-bcc triple points are artefacts of the approximation. Upon increasing 
the repulsion softness further and further, a stability gap eventually opens at T = for both fee 
and bee solids, which will be filled by one or more low-coordinated crystals. In the inset, the four 
potentials are plotted as continuous lines, with the respective rf(r) also shown as dotted lines. 

increasing pressures is (see js^] for acronyms) 

r 0.66 , 2.51 oa 3.50 , 3.86 5.07 oa 7.22 T1 „ 8.33 r 16.7 , 64.8 r 282 , 518 r 

ice — > bec — > pbn — > sn — > sc — y pbn — > cllo — > tec — > hep — > tec — > hep — > ice 

(4) 

where the numbers are the transition pressures in units of e/a 3 . We note that the internal 
parameter of the stable cI16 crystal varies from 0.043 to 0.047, incidentally close to the cI16 



phases of Li [38] and Na 39]. 



We have performed extensive Monte Carlo (MC) and molecular-dynamics simulations, 
using samples of 500 to 800 particles (simulations with 2048 particles lead to practically the 
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FIG. 2: (Color online). Theoretical melting lines (see text) for uyk{i")- Four cases are shown: 
a = 6.9 (blue), 5.7 (cyan), 4.5 (magenta), 3.3 (red). Repulsion softness grows in the arrow direction. 
Upon increasing the repulsion softness further and further, a stability gap eventually opens at T = 
for both fee and bec solids, which will be filled by one or more low-coordinated crystals. In the 
inset, the four potentials are plotted as continuous lines, with the respective rf(r) also shown as 
dotted lines. 

same results). Coexistence lines are computed by exact free-energy methods (details of the 



MC simulation are as in 25j). The phase diagram for P < 3 is plotted in Fig. 3. The bec 
melting temperature shows a maximum at P ~ 0.75, hence at not too low temperatures 
the bec solid melts upon compression into a denser fluid. At still lower temperatures, the 
bec phase transforms into a /3Sn phase. In the reentrant-fluid region a density anomaly 
occurs, i.e., the number density decreases upon cooling at constant pressure. This region 
is bounded from above by the temperature of the maximum density line (see Fig. 3). In 
Fig. 3, we also show points of diffusivity minimum and of — S2 maximum, S2 being the two- 



body entropy 



401 ] . Similarly to water, the region of structural anomaly encompasses that of 



anomalous diffusion, which in turn encloses the region of density anomaly. A compendium 
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FIG. 3: (Color online). Phase diagram of the YK potential for a = 3.3 {P and T are in reduced 
units). Full dots are two-phase coexistence points. Open dots are points of density maximum in 
the fluid phase. Diamonds and triangles denote points of —S2 maxima and D minima respectively 
(D being the self-diffusion coefficient), giving the left boundary of the regions of structural and 
diffusion anomaly (the right boundaries, defined by — S2 minima and D maxima, are out of the 
P range shown). Numerical errors are smaller than the symbols size. Within the pressure range 
2.5 -f- 3, the aspect ratio of the /3Sn phase is c/a ~ 0.60. Inset: uyk(^) for a = 3.3 (black solid 
line), f(r) = —u' YK (r) (blue dashed line), r/(r) (green dotted line), u'y K (r) (red dot-dashed line). 

of the anomalous behavior of the present system is shown in Fig. 4. 
III. DISCUSSION 

The anomalies exhibited by the YK fluid for a = 3.3 are quite similar to those of CS 
potentials. However, the much weaker softness of the former leads to peculiar local arrange- 
ment of particles in the fluid phase, as evidenced in the radial distribution function g(r). As 
P increases at constant temperature, the NN peak of g(r) moves gradually toward small r 
(Fig. 5, top panel). Meanwhile its height first grows, due to increasing proximity with the 
bcc solid, and then goes down in the pressure range where reentrant melting occurs. As 
P increases further, the NN peak of g(r) grows again while its position changes less and 
less sensibly due to the steep small-r repulsion. At low pressures, the behavior of g(r) is 
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FIG. 4: Thermodynamic and structural quantities for the YK fluid with a = 3.3. Left column: 
thermal expansion coefficient otp (units of kpjt), isothermal compressibility Kt (units of cr 3 /e), and 
constant-pressure specific heat Cp (units of kp) as a function of T along the isobar P = 2.5. For 
conventional liquids, otp, Kt, and Cp monotonically increase with T and sp > 0. Right column: 



translational order parameter — S2 (units of kp), bond-order parameter Qq [4l|], and self-diffusion 
coefficient D (units of a(e/m) 1 ^ 2 , where m is particle mass) as a function of P along the isotherm 
T = 0.06. For conventional liquids, — S2 and Qq increase with P while D decreases monotonically. 



similar to that of the Gaussian-core model 



42] (GCM), a fluid with a bounded interaction 



potential that also exhibits a bunch of waterlike anomalies. However, due to the absence of 
a true particle core, the behavior is completely different at high pressures where the GCM 
more and more resembles an "infinite-density ideal gas" with a radial distribution function 
g(r) = 1. To our knowledge, the g(r) behavior of the system investigated here has never 
been reported for other systems of truly impenetrable particles. 

To make clear the difference between the interaction studied here and the CS potentials 
investigated so far, we report the radial distribution function of a CS potential (Fig. 5, 
bottom panel). In this case, the heights of the first two peaks of g(r), associated respectively 
with the hard and soft length scales, change in opposite directions on increasing pressure, 
signalling that the "hard" NN distance becomes more and more populated at the expenses of 
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FIG. 5: Top panel: pair distribution function g(r) of the YK potential for a = 3.3. Bottom panel: 
g{r) of the YK potential for a = 2.1 (not far from the core-softening threshold a = 2.3). For 
a = 3.3, lines correspond to T = 0.06 and P = 0.1, 0.5, 1, 1.5, 2, 2.5, 3, 4, 6. For a = 2.1, curves refer 
to T = 0.07 and P = 0.5, 0.7, 1, 1.6, 2.3, 3, 4.1, 5.2, 6.4. The arrows mark the direction of pressure 
increase. 

the "soft" distance, while the position of the two peaks remains essentially unaltered. This 
behavior is a clear evidence of the simultaneous existence of two populations of particles 
having distinct effective diameters, a phenomenon that does not occur in the case analyzed 
here. 

To gain further insight we analyze the P dependence of the NN peak position r NN for 
the system under study, as compared with a CS potential {uyk(J' ) with a = 2.1) and u (r) 
(see Fig. 6). For the CS potential, r NN (P) consists of two branches, i.e., in a range of 

9 



NN 



1.2 - 



1.6 




0.8 







2 



4 



6 



P 



FIG. 6: Position tnn of the NN peak of g(r) in units of a as a function of P at constant T for: 
^yk(^), a = 2.1 and T = 0.07 (blue solid line and triangles); Uyk(0> a = ^>3 an< ^ = ^-06 (red 
solid line and full dots); uo(r), T = 0.06 (black dotted line, stopping near the melting point); uo(r), 
T = 1 (black dashed line). The gray scale is proportional to the height of the g(r) peak. The blue 
solid line without symbols represents a weighted average of the two r-NN branches, with weights 
proportional to the respective g(r) peak heights, for the case UYK.( r ) with a = 2.1. 

pressures two distinct populations of particles with different effective radii coexist, which 
legitimates the two-state-fluid picture. On the contrary, for the present system r NN (P) has 
only one branch (similarly to Uo(r) where, however, the g(r) peaks grow monotonously with 
P), meaning that all particles have the same effective radius. While the tnn(-P) curve for 
uo(r) is upward concave everywhere, in the present case it shows a portion with downward 
concavity, roughly corresponding to the reentrant-melting pressure range. The same feature 
is displayed by the weighted average of the two r NN branches for the CS potential (the 
downward-concave region is here shifted, as well as the reentrant portion of the melting line, 
towards smaller pressures). 

Thus, in the case considered here there exists only one type of local structure with an 
effective length scale that, however, shrinks with pressure in such a way as to capture the 
P-dependence of the mean length scale for the CS case. In particular, the shrinking rate 
shows a local maximum in the reentrant-melting pressure range, which makes the length 
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scale more loosely defined with respect to the adjacent regions. As shown by our results, 
this is sufficient to give origin to anomalous phase behavior. 



IV. CONCLUDING REMARKS 

Competition between two different local particle arrangements, arising from either direc- 
tional or core-softened isotropic forces, has so far deemed to be responsible for anomalous 
thermodynamic behavior. In this article, we have shown that the same behaviors may oc- 
cur also for weakly-softened (WS) potentials, i.e., simple fluids characterized by a repulsion 
that is only marginally softened and yields a single structure at a local level. Our results 
overturn the two-state-fluid picture derived from network-forming fluids and simplify the 
minimal scenario for the occurrence of anomalous behaviors, thus considerably widening the 
class of interactions that give rise to waterlike anomalies. WS potentials can be relevant 
in the realm of soft matter, where engineering interparticle forces is possible, and also for 
"hard" matter under extreme conditions. In the latter case, pressure typically induces suc- 
cessive rearrangements of the crystal structure in order to minimize the electronic energy. 



This is reflected at higher T in reentrant melting and (though more rarely reported |43j |) 
waterlike anomalies. The effect of pressure on the atomic structure may be dramatic, like 
e.g. in Cs where a 6s-5d band crossing causes a sudden drop in the effective radius of the 
atom 44j ("orbital collapse"), or less sharp as can be expected for lighter and smaller chem- 
ical elements. While two-scale CS potentials offer a simple model for systems undergoing 
orbital collapse, one-scale WS potentials may be appropriate to describe those elements 
where the atomic radius shrinks more gradually under pressure. 
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